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Abstract — The  power  spectrum  of  a  zero-mean  stationary  Gaus¬ 
sian  random  process  is  assumed  to  be  known  except  for  one  or 
more  parameters  which  are  to  be  estimated  from  an  observation 
of  the  process  during  a  finite  time  interval.  The  approximation  is 
introduced  that  the  coefficients  of  the  Fourier  series  expansion  of 
a  realization  of  long-time  duration  are  uncorrelated.  Based  on  this 
approximation  maximum  likelihood  estimates  are  derived  and 
fundamental  limits  on  the  variances  attainable  are  found  by  evalua¬ 
tion  of  the  Cram6r-Rao  lower  bound.  Parameters  specifically 
considered  are  amplitude,  center  frequency,  and  frequency  scale 
factor.  Also  considered  is  ripple  frequency  which  refers  to  the 
cosine  factor  in  the  spectrum  produced  by  the  addition  of  a  delayed 
replica  of  the  random  process.  The  dual  problem  of  estimating 
parameters  of  the  time-varying  power  level  of  a  nonstationary 
band-limited  white  noise  process  is  examined. 


I.  Introduction 

mHE  PROBLEM  of  estimating  parameters  of  the 


I  power  spectrum  of  a  stationary  random  process 
from  a  record  of  limited  duration  arises  in  a  variety 


of  applications.  For  example,  radio  astronomers  have 
mapped  the  structure  of  our  galaxy  by  measuring  the 
Doppler  shift  of  the  hydrogen  line  in  various  directions 
in  the  galactic  plane  [1].  In  radar  investigations  of  the 
ionosphere,  the  positive  ion  temperature  has  been  deduced 
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from  the  spectral  width  of  the  backscattered  signals  [2], 
Further  examples  pertaining  to  radar  astronomy  and 
seismology  are  described  in  Sections  VIII  and  X  of  this 
paper.  Some  of  the  methods  of  statistical  estimation 
theory  are  applied  here  to  these  problems.  With  the 
assumption  that  the  process  is  Gaussian  and  the  intro¬ 
duction  of  certain  approximations  valid  for  a  long  obser¬ 
vation  time,  an  expression  is  given  for  the  likelihood 
function  of  the  parameters.  From  this  expression,  maxi¬ 
mum  likelihood  estimates  can  be  determined  although 
they  cannot,  in  general,  be  stated  explicitly.  Funda¬ 
mental  limits  on  the  variance  attainable  by  these  or  any 
other  estimates  are  provided  by  evaluation  of  the  Cram6r- 
Rao  lower  bound.  Specific  results  are  obtained  for  the 
parameters  of  amplitude,  center  frequency,  frequency 
scale  factor,  and  ripple  frequency.  It  is  shown  that 
analogous  results  apply  to  the  estimation  of  parameters 
governing  the  time  variation  of  the  power  level  of  a 
nonstationary  band-limited  white  noise  process  such  as 
would  be  obtained  from  the  reflection  of  a  narrow  radar 
pulse  from  a  diffuse  cloud  of  scatterers. 

Many  authors  have  examined  the  problem  of  measuring 
the  overall  shape  of  a  power  spectrum  [3]-[5].  However, 
for  situations  like  the  above,  in  which  the  spectrum  is 
assumed  to  be  known  except  for  one  or  more  parameters, 
only  a  few  treatments  have  appeared  in  the  literature. 
Kelly,  Lyons,  and  Root  [6]  and  other  authors  analyzed 
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the  radiometer,  which  is  a  device  for  measurement  of 
spectral  amplitude.  Swerling  [7]  and  Mullen  [8]  described 
methods  for  estimating  frequency  scale  factor  and  center 
frequency  -which*  require  that  the  parameter  be  in  the 
vicinity  of  some  nominal  value  *and  that  the  amplitude 
of  the  power  spectrum  be  known  or  stabilized  by  an  auto¬ 
matic  gain  control.  They  obtained  expressions  for  the 
error  variances  which  for  power  spectra  symmetric  about 
a  center  frequency  are  equivalent  to  (37)  and  (41)  be¬ 
low.  Bogert,  Healy,  and  Tukey  [9]  discussed  a  method 
for  the  measurement  of  the  time  delay  of  an  “echo” 
of  the  random  process,  but  presented  no  theoretical 
evaluation. 

In  the  statistical  literature,  there  are  some  specialized 
results  dealing  with  parameters  of  moving  average  and 
autoregressive  discrete  time  series.  The  only  treatment 
which  applies  to  a  general  class  of  parameters  is  the 
work  of  Whittle  [10] — [12]  which  has  not  received  atten¬ 
tion  proportional  to  its  importance.  He  treated  hypothesis 
testing  for  a  discrete  Gaussian  time  series  with  unknown 
power  spectrum  parameters  by  approximating  the  co- 
variance  matrix  of  the  time  series  with  one  which  is 
circular  (periodic).  The  hypothesis-testing  framework 
places  the  amplitude  parameter  in  a  special  category 
and  leads  to  the  conclusion  that  its  maximum  likelihood 
estimate  is  asymptotically  uncorrelated  with  that  of  the 
other  parameters.  Here  (in  Section  IY)  a  different  result 
is  obtained  which  has  been  verified  by  other  means. 
Otherwise,  the  application  of  our  methods  to  discrete 
time  series  produces  results  consistent  with  Whittle’s. 

II.  Maximum  Likelihood  Estimates  and  the 
Cramer-Rao  Lower  Bound 

Some  general  results  of  statistical  estimation  theory  [13] 
are  presented  in  this  section.  Consider  an  iV-dimensional 
vector  random  variable  £  with  probability  density  g(x;  a) 
depending  on  a  parameter  a  with  true  value  a0.  For  an 
observation  of  £,  the  natural  logarithm  of  the  likelihood 
of  a  is  defined  as 

A(£;  a)  =  log  g&;  a) 

The  maximum  likelihood  (ML)  estimate  &  is  the  value  of 
a  which  maximizes  A (5;  a).  It  can  sometimes  be  found 
explicitly  as  the  root  of  the  likelihood  equation 


(i) 


but  usually  a  linearization  in  the  vicinity  of  a0  or  a 
method  of  successive  approximations  is  required. 

The  Cram6r-Rao  lower  bound  states  that  under  general 
regularity  conditions  for  any  estimate  a, 

[1  +  db/da0]2 

-ElC  (?;a0) 


vara,  > 


where 


b  =  Ea,  —  aQ 


is  the  bias  of  a,.  An  alternate  expression  for  the  denomi¬ 
nator  of  (2)  is  obtained  from  the  identity 

-  EA''(Z;  a0)  =  £[A'(f;  «„)f  (3) 

In  this  paper,  a  prime  always  denotes  differentiation  with 
respect  to  a  and 


A'(5;  a0)  =  ^  A(?;  a) 


A"«;  «o)  =  ^5  Afo  a) 


-MU 


An  unbiased  estimate  whose  variance  satisfies  (2)  with 
the  equals  sign  is  said  to  be  efficient.  However,  efficient 
estimates  exist  only  in  certain  cases.  ML  estimates  are 
not  necessarily  unbiased  or  efficient.  However,  ML  esti¬ 
mates  based  on  N  independent  samples  of  a  random 
variable  have  certain  optimal  properties  as  iY  — *  <». 
Grenander  [14]  extended  these  results  to  ML  estimates 
based  on  a  realization  of  a  random  process  of  duration 
T.  He  established  under  general  conditions  that  as  T  — >  oo  f 
&  converges  in  probability  to  an  and  is  asymptotically 
efficient.  This  implies  that  V  —  EA"fe;  a0)(d  —  oco)  con¬ 
verges  in  distribution  to  a  distribution  with  mean  zero 
and  variance  one.  In  addition,  this  limiting  distribution 
is  usually  Gaussian. 

Frequently  the  bias  of  an  estimate  under  considera¬ 
tion  cannot  be  established,  so  (2)  is  not  informative. 
However,  the  Cram6r-Rao  lower  bound  can  still  be  in¬ 
terpreted  in  terms  of  the  “sensitivity”  introduced  by 
Kelly,  Lyons,  and  Root  [6].  They  consider  any  statistic 
a4  which  is  a  measure  of  aQ  in  the  sense  that  its  expecta¬ 
tion  is  a  monotonic  function  of  a0  and  define 


sensitivity  = 


dEaJdot  o 


standard  deviation  of  a,  for  a  =  a0 


The  reciprocal  of  the  sensitivity  is  just  that  small  change 
in  a0  required  to  change  the  mean  value  of  a,  by  one 
standard  deviation.  It  is  seen  that  (2)  can  now  be  written 


G 


1 


1 


ensitivity/ 


EA  "(?;a0) 


(4) 


Consider  now  the  general  class  of  estimates  a *  ob-  1 
tained  by  maximizing  over  a  some  function  T((;  a) 
which  depends  on  £  and  a  (but  not,  of  course,  on  «0) 
and  has  the  further  property 


Er'£;a0)  =  0 


(5) 


The  dispersion  of  the  sampling  errors  can  be  character¬ 
ized  by  the  quantity 


(2) 


«,  *  E[r'&ao)Y 
Sia  *  - 
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This  quantity  can  be  further  interpreted  when  r(£;  a)  is 
sufficiently  regular  so  that  in  the  vicinity  of  the  true 
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parameter  value  a0  it  can  be  approximated  by  the  first 
5jv  three  terms  of  the  Taylor  series  expansion 
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rfo«„)  +  (a  -a0)r'($;a0) 
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(a  ~  «o)2 


r"(?;  «o)  (6) 


tat'' 

The  maximum  value  of  I’(£;  a)  over  a  occurs  where 

li^k.  =  o 

.  so  from  (6) 
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S(a*) 
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(10) 


The  analysis  is  based  on  the  following  assumptions: 


1)  The  random  process  is  stationary,  Gaussian,  and 
zero  mean  with  a  double-sided  power  spectrum  P(j)  a ) 
which  is  a  known  function  of  the  parameter  a  whose 
value  is  to  be  estimated.  (If  the  process  consists  of  a 
random  signal  plus  an  independent  noise,  then  P(j\  a)  is 
sum  of  the  two  individual  spectra.) 

A  narfinnkr  r^nlirzatirm  v(t\  nf  t 


<  a0  <  am 


outside  of  ft 

fcwc  *-•  1 


a 

Sisk 


a  finite-dimensional  formulation  of  the  problem  and  is 
meaningful  in  the  usual  practical  situation  where  fre¬ 
quencies  above  some  finite  limit  are  not  observable.  The 
assumption  fx  >  0  is  only  for  the  convenience  of  elimi¬ 
nating  zero  frequency  which  enters  unsymmetrically 
into  the  analysis.) 

5)  For  all  /  and  all  a  such  that  amin  <  a  <  am«  the  func¬ 
tions  P(j;  a),  dP(j;  a)/dj,  dP(f ;  a)/da,  and  d2P(j ;  a)/df  da 
are  bounded  and  continuous.  Also  P(/;  a)  is  bounded  away 
from  zero  for  fx  -  (/2  -  fx)  <  /  <  /2  +  (/2  —  /i). 


(8) 


Then  the  distribution  of  the  sampling  error,  a*  —  a0, 
is  approximately  the  same  as  that  of  the  expression  (8). 
If  the  further  approximation  is  made  that  T"({;  a0)  is  a 
constant  or  that  its  random  components  are  of  a  smaller 
order  than  those  of  the  numerator  of  (8),  then  the  random 
fluctuations  of  r'(£;  a0)  about  zero  may  be  thought  of 
as  a  linear  noise  term  displacing  the  peak  of  the  parabolic 
approximation  to  r(£;  a)  and 


In  dealing  with  statistical  inference  for  random  proc¬ 
esses,  a  basic  problem  is  to  choose  a  set  of  “observable 
coordinates”1  which  represent  the  process  and  whose 
probability  distribution  is  tractable.  In  this  paper,  the 
observable  coordinates  are  taken  as  the  Fourier  coefficients 


=  _i_  f 

Vf  Jo 


v{t)e~ 


dt 


(9) 


It  has  been  shown  by  Godambe  [15]  that  if  S(a*)  is 
taken  as  a  measure  of  the  dispersion  of  the  estimate  a*, 
without  necessarily  referring  to  the  particular  interpre¬ 
tation  mentioned  above,  then  under  general  regularity 
conditions  £(<**)  is  minimized  when  T(^;  a)  =  A({;  a) 
and  a*  =  d.  This  minimum  value  is,  with  (3), 

m  .  *[v«;->r  _ 


where  /0  =  1/T.  This  is  plausible  since  under  general 
conditions,  for  almost  every  realization,  v(t)  can  be  repre¬ 
sented  by  the  limit  in  the  mean  of  its  Fourier  series 
expansion  [16].  The  real  and  imaginary  parts  of  the  yn 
have  a  multivariate  Gaussian  distribution  and  the  likeli¬ 
hood  function  depends  on  the  covariance  matrix  of  the 
entire  set  of  coefficients  which  is  too  complicated  to  be 
useful.  We  now  observe  [16],  [17]  that  with  the  above 
definition 


lim  Eynyn  =  F{n/T) 


n  =  m 


(id 


(0 


n  m 


This  same  quantity  appears  in  the  Cram6r-Rao  lower 
bound  (2),  but  the  present  result  holds  with  the  equals 
sign  for  any  regular  estimate,  biased  or  not. 


where  n  — >  °o  as  T  — >  «>  so  that  n/T  remains  constant 
and  ym  is  the  complex  conjugate  of  ym.  It  can  be  shown 
wTith  the  above  assumptions  that  this  convergence  goes 
uniformly  as  log  T/T.  The  approximation  is  now  intro¬ 
duced  that  for  large  but  finite  T, 


III.  Approximate  Maximum  Likelihood  Estimates  for 
Power  Spectrum  Parameters 


Py  nY  m  — 


n  =  m 


(12) 


n  m 


2)  A  particular  realization  v(t)  of  the  process  is  ob¬ 
served  for  0  <  t  <  T. 


3)  The  true  value  of  a  is  denoted  by  a0  and  the  actual 
power  spectrum  is  <t>(/)  =  P(/;  a0).  a0  is  assumed  to  lie 


‘  within  a  known  finite  interval  an 

4)  For  any  a  within  this  interval,  only  values  of  P(j;  a) 
within  some  fixed  finite  range  0  <  /i  <  |/|  <  <®  depend 
v  on  a.  Since  fx  and  /2  may  be  anywhere  outside  of  the 
range  of  dependence,  for  convenience  we  take  fx  =  Nx/T 
and  /2  =  N2/T  where  Nx  and  N2  are  integers.  (This 
assumption  means  that  P(j;  a)  has  some  constant  shape 
<  |/|  <  /2.  The  assumption  /2  <  00  allows 


The  Karhunen-Lofcve  expansion2  which  is  frequently 
employed  in  the  analysis  of  random  processes  is  not  well 
adapted  to  spectral  parameter  estimation  problems  since 
it  is  in  terms  of  an  uncorrelated  set  of  eigenfunctions 
which,  except  in  special  cases,  change  in  a  complicated 
way  as  the  parameters  vary.  However,  if  and  only  if  the 
covariance  function  of  the  random  process  is  periodic 
with  period  T  then  the  eigenfunctions  become  complex 
exponentials,  the  eigenvalues  are  values  of  the  power 
spectrum  spaced  at  intervals  l/T,  and  the  coefficients  of 
the  Karhunen-Loeve  expansion  are  uncorrelated'  [16]. 
Therefore,  the  approximation  (12)  is  equivalent  to  ap¬ 
proximating  the  covariance  function  of  the  process  by 
one  which  is  extended  periodically. 

The  likelihood  function  of  the  parameter  a  is  now 


1  Grenander  [14],  pp.  207-209. 

2  pp.  199-200;  Davenport  and 
96-101. 


Root,  op.  cit.  [17],  pp. 
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for  evaluating  the  Cram6r-Rao  lower  bound.  In  fact  (17) 
and  (18)  converge  at  least  as  rapidly  as  (log  T)4/Ti 
Experiences  of  other  investigators  suggest  that  good 

results  are  obtained  for  modest  values  of  T.  but  a  nu- 

7 

merical  guide  is  not  yet  available. 

The  step  leading  from  (8)  to  (9)  can  also  be  justified 
so  this,  together  with  (10),  (16),  (17),  and  (18),  enables 
us  to  state  the  basic  conclusion  that  for  large  T 


1965 

obtained  by  substituting  the  observed  values  of  the 
Fourier  coefficients  into  the  complex  Gaussian  prob¬ 
ability  density  [18]  which  with  (12),  is 


.V: 


1 


7T P(n/0;  a) 


exp 


Jin 


P(nfQ;  a)_ 


and  depends  upon  v{t)  only  through  the  “periodogram” 
|7n|2.  The  approximate  ML  estimate  a *  is  found  by 
maximizing  the  approximate  log  likelihood 


rA 


w 


r 

r(r;«)  =  -  12 

»-jv i  L 


log  7rP(n/0;  a)  + 


!?» 

P(nj o 


(13) 


r(r; «) 


log  irP(j;  a )  + 


P(i; 


dj  (14) 


r'(Y-  a)  =  -  y. 

{y,)  rw;a)  i 


P2(nfa 


-ET"(r,«o) 


'(n/0;  a0)  p 
_P(nf0)  a0)  J 

['-  rP'(/;ao)T 
I.  Lp(/;  «o)  J 


df 


i ; ,  L  S(«*) 

lini  n/ a\ 

r-.»  <S(<S) 

It  can  also  be  shown  that 


=  1 


J.™  EA"(y,  a.) 


var  a’ 


If  P(/;  a)  is  a  slowly  varying  function,  then  the  summa¬ 
tion  (13)  can  be  approximated  by  an  integral.  The  indi¬ 
vidual  |7n|2  need  not  be  measured  since  their  values  can 
be  smoothed  and  replaced  in  this  integral  by  a  continuous 
function  <£*(/)  which  is  equivalent  to  an  approximately 
unbiased  overall  spectral  estimate  such  as  given  by  the 
lackman-Tukey  method  [4]  or  by  a  power  spectrum 
analyzer.  For  the  integral  to  be  a  good  approximation, 
all  the  |7n|2  must  be  effectively  included  and  the  resolution 
must  be  sufficient  so  that  the  detailed  structure  of  the 
spectrum  is  not  obscured.  Then 


(15) 


(16) 


~  l M  l_P(n/o; «.)  J  f 

“{'rtw&N" 

-wrii106™ 


;«)T  djV 

Ja-dt  ) 
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and  these  expressions  are  also  asymptotic  values  for  the 


Cram6r-Rao  lower  bound. 


IV.  Joint  Estimates 


The  accuracy  of  these  estimates  can  be  characterized 
in  terms  of  the  Cram6r-Rao  lower  bound  by  approxi¬ 
mating  PA"(y;  «o)  with  ET"(y;  a0).  We  have 

|2p, 


For  a  set  of  M  unknown  parameters  a{,  denoted  by 
the  vector  a,  joint  ML  estimates  are  found  by  maximizing 
A(£;  a)  over  a.  There  is  a  corresponding  extension  of  the 
Cram6r-Rao  lower  bound  for  any  set  of  joint  unbiased 
estimates  a,.  The  covariance  matrix,  denoted  by  [cov  (a.)], 
has  elements 


COV  oteijati  =  E[aei  —  Eati][ati  —  Ea0j) 
Let  [<£]  be  the  matrix  with  elements 


where  P'(nf0',a)  =  ( d/da)P(nj0;a ).  It  is  noted  that 
£r'(y;a0)  =  0,  satisfying  (5).  A  second  differentiation 
and  the  substitution  of  (3)  give 


£if  =  -  E 


A(£;  «) 


datd<Xj 

Then  the  matrix 

[cov  («.)]  -  [£]-’ 

is  positive  semidefinite.  This  implies  that 
var  >  l/£u 


9 

•■Y  S-lH! 


(20) 


How  good  is  this  approximate  analysis?  Although 
similar  approximations  have  been  employed  by  Bryn  [19], 
Good  [20],  Freiberger  [21],  and  others,  as  well  as  by 
Whittle,  a  quantitative  justification  has  not  been  given. 
We  have  been  able  to  establish  by  somewhat  involved 
arguments  [22]  (whose  length  discourages  inclusion  here) 
that  with  assumptions  l)-5)  above,  the  approximate  ML 
estimates  are  as  good  as  the  exact  ML  estimates  in  the 
sense  that 


for  i  —  1,2,  •  •  •  ,  M.  Analogous  results  can  be  work 
out  for  the  sensitivity  and  the  spread. 

For  the  power  spectrum  parameter  estimation  prob¬ 
lem  we  find,  corresponding  to  (16), 


/ 


-rf 

J  r  i 


A  dP(J ;  a)  dP(j ;  a) 


dcti 


da* 


dj 


(21) 


Sf.  J'fi) 


(17) 


(18) 


A  case  which  is  important  in  practical  applications  is 
that  in  which  the  amplitude  level  A  of  the  spectrum  is 
one  of  the  unknown  parameters,  so  that  the  spectrum 
has  the  form  AP(f;  a).  Frequently  A  is  a  nuisance  param¬ 
eter  whose  value  is  not  required.  This  often  occurs  in 
electronic  equipment  where  the  absolute  value  of  the 
gain  is  difficult  to  calibrate.  For  a  fixed  value  of  a  (see 
Section  V), 


so  (16)  furnishes  a  convenient  asymptotic  approximation 


A*  = 


1 


yi 


T(U  -  /,)  „  P(njo  J  a) 


- - 

7? *  .‘i  z'\i  .  n 
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Substituting  this  expression  back  into  r(A,  a)  gives  the  This  is  an  unbiased  estimate  and  the  variance  is  found 
log  likelihood  in  terms  of  a  alone.  From  this  it  can  be  to  be  exactly 
shown  that  a *  can  be  obtained  by  maximizing  over  a 


r(x4*,  a)  =  —  log  xP(n/0;  a) 

h  T(l,  -  /,)  log  [?P(i);j”«)]} 

-  il, 


var  A* 

Al 


1 


T(J ,  -  /,) 


(28) 


t\£,  ■ 

-  ^  V  . 

F?  1 


+  constant 


log:rP(/;a)d/ 


Therefore,  as  (/2  —  /0  approaches  infinity  the  variance 
approaches  zero  for  any  fixed  T .  This  result  would  be 
expected  from  the  well-known  singularity  of  this  situa¬ 
tion  [23]. 

A  second  special  case  occurs  when  AS(J )  <3C  N0  for 
all  /.  A*  can  be  approximated  by  setting 


? ,  .  i  \ 


+ 


(U  -  /,)  log  [r  J,m  df  }  +  constant  AS(nf0)  +  N0 

(22)  in  the  denominator  of  (25)  so 


No 


without  the  necessity  of  determining  A *. 

The  effect  of  the  lack  of  knowledge  of  A  is  placed  in 
evidence  by  evaluating  the  approximation  to  the  co- 


r'(r;  A) 


Sjnfo) 

Nl 


[ 


A  S(nj0)  +  N0-  |T, 


(29) 


variance 
to  var  a*  is 


matrix  of  a*  and  A*.  The  element  corresponding  SettinS  (29>  eclual  to  zero  8ives 


var 


SV- 


-“KM1* 

j  -mh>[C»lT  « 

When  this  is  compared  with  (19)  it  is  seen  that  the  lack 
of  knowledge  of  A  increases  var  a*  unless  the  second 
integral  in  (23)  is  zero. 

v 

V.  Estimation  of  Spectral  Amplitude 
The  power  spectrum  is  assumed  to  be 


r  [**(/)  -  Km  df 

i*  ~  ih _ 


f 


(30) 


S\f)  df 

This  estimate  is  also  unbiased  and 
var  A* _ Nl 


Al 


T  f''  [^0^(/)]a  df 


(31) 


The  estimates  (27)  and  (30)  are  both  linear  functionals 
of  $*(/)  and  can  be  obtained  as  linear  functions  of  the 
outputs  of  suitable  radiometers.  A  radiometer  consists  of 
a  filter  with  transfer  function  //(/)  followed  by  a  square 
law  detector  and  integrator.  To  evaluate  the  estimate 
where  <S(/)  and  N0  are  known  and  the  parameter  A  is  to  (27)  it  is  required  that  |//(/)|2  =  1  /S(f)  while  for  (30), 


P(f;  A)  =  AS(J)  +  No 


be  estimated.  Then  from  (13) 


pV, 


r(r;  A) 

i 

l 


N>  f 

=  —  S  log  71- 

n-Af,  V 


[A  S(nf0)  +  N0] 

hi 


+ 


MS(n/„)  +  N 


and 


1,  n  v  S(nfo)\AS(nf0)  +  N0  -  |?J2] 

r  <T’ A)  •  ? _  hswj  +  >/.r 


(24) 


(25) 


var  A  * 


From  (19) 


|//(/)|2  =  £(/)•  The  variances  (28)  and  (31)  agree  with 
known  radiometer  formulas. 

VI.  Estimation  of  Center  Frequency 
In  this  case,  the  power  spectrum  has  the  form 

P(/;/<)  =P(/-/,)  />o 

and 

r(r;  /.)  =  -  L  log  rP(nf0  -  fe)  -  j]  p  \7^  — 

n-Ni  n-Nt  P \Jl] Q  J c) 


1 


r  Sink)  I2!- 

«  -  T  ^  log  irPif-fJdf  -T  f 

Uo<S(n/0)  +  N0]  J 

"  f  i  "  ) 

^v" 

w:  :  !  ; 


L  [’ 

■  r  sin  t 

\TL 

lA0S(f)  +  No] 

(26) 


Let  us  further  assume  in  this  section  that 

$(/i)  =  $(/2)  =  $0 


(33) 


A  general  explicit  solution  for  A*  is  not  known.  How-  for  all  allowable  /c  so  that  $(/)  has  the  form  indicated  in 


ever,  there  are  two  significant  special  cases.  The  first  is  Fig.  1.  Then  the  first  integral  in  (32)  is  independent  of 


for  N0  =  0.  Then  setting  r'(Y;  A)  =  0,  we  obtain 

i  rfm  **v) 


fc  and  fe*  is  obtained  by  minimizing  over  fe 


-j—  r 

U  -  U  in 


S(f) 


df 


(27) 


in  P(f  - 


/.  P(f  ~  /<) 


df 


(34) 
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Fig.  1.  The  assumed  spectral  shape  for  the  estimation  of  Doppler 

shift. 
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f 


**(/)?(/  -  /«)  df 
Var  f*  is  now  evaluated  from  (19).  For  /  >  0, 

P’V  -/<)  =  - 

and  (19)  becomes 


VII.  Estimation  of  Frequency  Scale  Factor 
(Spectral  Width) 

The  power  spectrum  is  assumed  to  be 


P(f;  h)  =  P[h(f  -  /,)]  f  >  0 


r(r;  h) 


-rlC 


log  -/,)]<!/ 

+  /' 


**(/)<*/  1 
W  -  /.)]J 


=  -  T[It  +  /,) 


(38) 


Fig.  2.  The  assumed  spectral  shape  for  the  estimation  of  frequency 

*  ! fact 


This  can  be  thought  of  as  a  “cross-correlation”  between 
<£*(/)  and  the  reciprocal  of  P(j  —  /«.).  Knowledge  of  the 
amplitude  level  of  the  spectrum  is  not  required. 

On  the  other  hand  if 


(35) 


where  |p(/  —  fe)\  <C  $0  for  all  /,  then  j*  is  found  approxi¬ 
mately  by  maximizing 


r(r;  h) 


(36) 


where 


K 


(37) 


The  features  of  the  spectrum  which  make  possible  the 
estimation  of  fe  are  measured  in  a  sense  by  the  integral 
in  (37).  The  integrand  is  non-negative  at  all  frequencies 
and  is  large  where  <£(/)  is  changing  rapidly.  From  (23)  it 
is  found  that  when  (33)  holds,  the  lack  of  knowledge  of 
the  spectral  amplitude  does  not  increase  the  approxi¬ 
mation  to  var  f*. 


var 


Now  assume  (33)  holds  for  all  allowable  h  so  that  $(/) 
has  the  form  indicated  in  Fig.  2.  Then  the  evaluation  of 
7i  can  be  simplified  by  writing 

P[h(f  ~  /«)]  =  *.{1  +  p[h(f  -  /,)]}  (39) 


By  a  simple  transformation  of  variables 


-'[M 


+  vii  -  U))  df 


The  maximization  of  (40)  to  determine  h*  is  seen  to 
depend  upon  the  knowledge  of  the  spectral  amplitude. 
If  it  is  unknown  then  the  method  of  Section  IV  can  be 
used. 

The  Cramdr-Rao  lower  bound  (19)  is  found  to  be 


**  *  (r  /: « -  dt/ «®j  «r  3  #> 


In  comparison  with  (37),  the  integrand  in  (41)  has  an 
additional  factor  (/  —  fe)2  which  emphasizes  the  spectral 
variations  by  the  square  of  their  distances  from  fe.  Also 
in  contrast  to  the  estimation  of  /c,  an  evaluation  of  (23) 
for  the  present  case  shows  that  when  the  spectral  ampli¬ 
tude  is  unknown,  var  h*  is,  in  general,  increased. 

r  -g* 

VIII.  Estimation  of  Ripple  Frequency  l 


A  problem  which  arises  in  seismology  and  other 
[9]  is  that  for  which 


where  fe  is  a  known  center  frequency  and  P[h(j  —  /«.)]  is 
known  except  for  the  frequency  scale  factor  h  which  is 
to  be  estimated.  (Note  that  the  Doppler  effect  actually 
is  a  change  of  scale  factor  with  fe  =  0  and  the  Doppler 
shift  is  an  approximation  valid  when  the  width  of  the 
power  spectrum  is  small  compared  with  its  center  fre¬ 
quency.)  h*  is  obtained  by  maximizing 


v{t)  =  A[u{t)  +  au(t  -  t)] 

where  u(t)  is  a  zero-mean  stationary7  Gaussian  random 
process  with  knowrn  power  spectrum  £(/),  A  and  a  are 
unknown  amplitude  factors,  and  r  is  a  time  delay  to  be 
estimated.  Then 


P(/;  A,  a,  r)  =  A$(/)[l  +  2a  cos  +  a2 J 

and  the  addition  of  the  “echo”  a  u{t  —  r)  produces  a 
cosine  ripple  of  frequency  r  in  the  power  spectrum.  It 
can  be  shown  that  r*  is  obtained  byr  minimizing  over 
r  and  a 


r(r;  A 


**(/) _ 

S(f)[]  +  2a  cos  2 irfr  +  a2]  ’  (  i) 


rt-.. 


r*v  i 
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For  small  a  this  is  equivalent  to  maximizing  over  r 

,  f'-ig£cm2rlTdl 
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2(TF2  -  IF,) 
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varz(nA/)  =  R(nAt;  0) 


(45) 
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This  has  a  form  almost  identical  with  (13)  and  the  maxi¬ 
mum  likelihood  estimation  of  0  is  analogous  to  the  dis¬ 
cussion  of  the  preceding  sections.  The  result  corresponding 
to  (19)  is 


which  amounts  to  “whitening”  <£*(/)  and  finding  the 
peak  of  the  cosine  transform.  For  |a|  <<C  1  and  /2  ^  l/|r| 
'="  it  is  found  that 

&  \?r  a*T3Ul  -  It)}  (44) 


var 


IX.  A  Dual  Problem 

+  i 

A  problem  which  is  dual  in  a  certain  sense  to  that 
stated  in  Section  III  is  now  examined.  Consider  a  zero- 
mean  nonstationary  Gaussian  random  process  y(jt)  with 
a  time-varying  covariance  function 

E[y{t)y{t  +  r)]  =  R(t;  fi)S(r) 

he.,  y(t)  is  a  white  noise  process  with  a  power  level  R(t\  0) 
which  is  a  function  of  time  and  of  a  parameter  0 .  Suppose 
R (t;  0)  is  known  except  for  0  which  is  to  be  estimated. 
The  following  assumptions  on  R(t;  0)  are  analogous  to 
those  of  Section  III  regarding  P(J ;  a)  with  time  and  fre¬ 
quency  interchanged:  Only  a  finite  frequency  range  of 
the  y(t)  process  can  be  observed  so  that,  in  effect,  it  has 
been  passed  through  an  ideal  band-pass  filter  with  re¬ 
sponse  limited  to  frequencies  Wx  <  f  <  W2  and  the 
estimates  must  be  based  on  the  resulting  band-limited 
output  process  yB(t).  For  any  value  of  0,  only  values  of 
R(t;  0)  within  a  finite  interval  tx  <  t  <  t2  depend  on  0. 
ti  and  t2  do  not  depend  on  0  and,  for  convenience,  we  take 


where  M x  and  M2  are  integers.  This  assumption  is  valid 
if  R(t;  0)  has  a  constant  shape  outside  of  or  is  only  observ¬ 
able  within  some  fixed  time  interval.  Finally,  R(t;  0) 
obeys  certain  general  regularity  conditions. 

Introduced  the  normalization 


Z(t)  =  yB(t)  Vw2  -  IF, 

On  the  basis  of  the  above  assumptions  the  band-limited 
process  z(t)  can  be  represented  by  the 

2(<j  -  f,)(IF3  -  IF,)  +  1 

samples  z(nAt)  where  At  =  hOV2  —  Wx)  and  n  =  MX) 
Mx  +  1,  •  •  •  ,  M2 .  The  z(nAt)  are  Gaussian  and  to  a  good 
approximation  are  independent  with 


Thus  the  samples  z(nAt)  are  duals  of  the  real  and  imagi¬ 
nary  parts  of  the  Fourier  coefficients  of  Section  III. 
The  logarithm  of  the  likelihood  of  the  set  of  samples  is 
approximately 

r(Z! »  -  - 1 .5.  [los  2'S(n4,: »  +  bS]  (46) 


(47) 

Thus,  the  wider  the  observable  bandwidth  (TF2  —  Wx), 
the  smaller  the  variance  of  the  estimate. 

X.  An  Application  to  Radar  Astronomy 

An  interesting  example  of  a  possible  application  of 
these  results  is  in  radar  astronomy  [24],  [25].  When  a 
CW  signal  is  transmitted,  the  energy  reflected  from  the 
moon  or  a  planet  can  be  considered  as  a  Gaussian  random 
process.  The  radial  velocity  of  the  body  can  be  inferred 
from  the  Doppler  shift  and  the  rotation  rate  from  the 
spectral  width  (frequency  scale  factor),  so  the  analyses 
of  Sections  VI  and  VII  apply.  The  effect  of  thermal 
noise  at  the  receiver  is  included  by  adding  its  constant 
level  to  the  power  spectrum  of  the  received  energy. 

On  the  other  hand,  Section  IX  applies  to  the  case 
where  a  narrow  pulse  of  bandwidth  ( W2  —  Wx)  is  re¬ 
flected  from  a  target  having  a  large  delay  spread.  If  the 
target  consists  of  many  point  scatterers  randomly  and 
independently  distributed  with  a  spatial  density  pro¬ 
portional  to  wR{t)  0 ),  then  the  received  signal  approxi¬ 
mates  a  Gaussian  nonstationary  band-limited  white 
noise  process  with  time-varying  power  level  R(t\  0).  The 
target  characteristics  can  be  estimated  as  the  parameters 
of  R(t;  0)y  the  range  being  analogous  to  the  center  fre¬ 
quency  and  the  target  depth  to  the  frequency  scale 
factor  which  were  discussed  in  Sections  VI  and  VII. 
Another  parameter  of  significance  is  the  rapidity  with 
which  the  power  level  falls  off  with  range,  since  it  provides 
a  measure  of  the  surface  roughness.  It  should  be  noted 
that  the  regularity  assumptions  may  not  be  satisfied  in 
some  cases  because  of  the  sharp  leading  edge  of  the 
return  signal. 
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